Optimal Jacobian accumulation is NP-complete

We show that the problem of accumulating Jacobian matrices by using a minimal number of floating-point operations is NP-complete by reduction from Ensemble Computation. The proof makes use of the fact that, deviating from the state-of-the-art assumption, algebraic dependences can exist between the local partial derivatives. It follows immediately that the same problem for directional derivatives, adjoints, and higher derivatives is NP-complete, too. We consider the automatic differentiation (AD) [2] of an implementation of a non-linear vector function y = F (x, a), F : IR ⊇ D → IR, as a computer program. With the Jacobian matrix F ′ of F defined as usual tangent-linear (ẏ = F (x, a) ∗ ẋ, ẋ ∈ IR) and adjoint (x̄ = (F (x, a)) ∗ ȳ, ȳ ∈ IR) versions of numerical simulation programs with potentially complicated intraand interprocedural flow of control can be generated automatically by AD tools. This technique has been proved extremely useful in the context of numerous applications of computational science and engineering requiring numerical methods that are based on derivative information. For the purpose of this paper we may assume trivial flow of control in the form of a straight-line program. Similarly, one may consider the evaluation of an arbitrary function at a given point to fix the flow of control. Our interest lies in the computation of the Jacobian of the active outputs (or dependent variables) y = (yj)j=1,...,m with respect to the active inputs (or independent variables) x = (xi)i=1,...,n. The ñ−vector a contains all passive inputs. Conceptually, AD decomposes the program into a sequence of scalar assignments vj = φj(vi)i≺j for j = 1, . . . , p + m, where we follow the notation in [2]. We refer to this equation as the code list of F, and we set xi = vi−n for i = 1, . . . , n and yj = vp+j for j = 1, . . . , m. The vj , j = 1, . . . , p, are referred to as intermediate variables. The notation i ≺ j marks a direct dependence of vj on vi meaning that vi is an argument of the elemental function 3 φj . The code list induces a directed acyclic graph G = (V, E) such that V = {1−n, . . . , p+m} and (i, j) ∈ E ⇔ i ≺ j. Assuming that all elemental functions are continuously differentiable at their respective arguments all local partial derivatives can be computed by a single evaluation of the linearized code list cj,i = ∂φj ∂vi (vk)k≺j ∀i ≺ j vj = φj(vi)i≺j Software and Tools for Computational Engineering, Department of Computer Science, RWTH Aachen University, 52056 Aachen, Germany, http://www.stce.rwth-aachen.de, [email protected] 2 F is used to refer to the given implementation. Elemental functions are the arithmetic operators and intrinsic functions provided by the programming language. for j = 1, . . . , p + m and for given values of x and a. The corresponding linearized version of G is obtained by attaching the cj,i to the corresponding edges (i, j). Various elimination techniques have been proposed for efficient Jacobian accumulation based on G [3]. Theorem 1 Optimal Jacobian Accumulation is NP-complete. We reduce from Ensemble Computation [1]. A given solution is verified in polynomial time by counting the number of operations. Given an arbitrary instance of Ensemble Computation we define the corresponding Optimal Jacobian Accumulation problem as follows: Consider y = F (x, a) where x ∈ IR, a ≡ (aj)j=1,...,|A| ∈ IR |A| is a vector containing all elements of A, and F : IR → IR defined as